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■ ABSTRACT 

Q ' We numerically investigate the dynamics of how a close stellar fly-by encounter of a 

symmetrical circumstellar planetesimal disc can give rise to the many kinds of asym- 
metries and substructures attributed to the edge-on dusty disc of (3 Pic. In addition we 
present new optical coronagraphic observations of the outer parts of (3 Pic's disc, and 
report that the radial extent is significantly greater than was found in previous mea- 
surements. The northeasterly extension of the disc's midplane is now measured out to 
1835 au from the star; the southwesterly component is measured out to 1450 au. Thus 
the asymmetry in the length of the former with respect to the latter is approximately 
| 25 per cent. We proceed to use the length asymmetry induced in the distribution of 

, simulation test particles as the principal diagnostic feature when modelling the disc 

response in order to constrain fly-by parameters. In particular we favour a low in- 
clination prograde and near-parabolic orbit perturber of mass approximately 0.5 M©. 
. These initial conditions suggest that the perturber could have been physically associ- 

ated with [3 Pic prior to the encounter. Thus we also consider the possibility that the 
r-j ! perturber could be bound to (3 Pic: a consideration also of general interest where dust 

O discs are known to exist in binary star systems. We show that a further consequence 

of a low velocity encounter is that the perturber could have captured planetesimals 
from the f3 Pic disc, and we deduce that as a result of this the perturber could dis- 
play a dust disc that is presently amenable to observation. In some of our models, 
we can relate groupings of perturbed particles to the large-scale structure of the (3 
Pic disc. The groupings correspond to: high eccentricity and inclination particles that 
reach apocentre and maximum height in the southwest, moderately eccentric and low 
inclination particles that reach apocentre in the northeast, and relatively unperturbed 
particles inside ^200 au radius. 
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1 INTRODUCTION 

The starlight of j3 Pic scatters off dust grains in its edge- 
on circumstellar disc (Smith & Terrile 1984). Optical imag- 
ing reveals that the disc exhibits asymmetrical structure 
with respect an ideal disc on projected radial scales from 
~50au (Burrows et al. 1995, Mouillet et al. 1997) to 
~1000au (Smith & Terrile 1987, Kalas & Jewitt 1995). 
Kalas & Jewitt (1995) presented sensitive _R-band images 
that showed the northeast (NE) extension could be detected 
out to 790 au, whereas the southwest (SW) extension was 
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detectable out to just 650 au, corresponding to a ~20 per 
cent length asymmetry. More recently, the longer NE mid- 
plane extension has been discovered to carry several bright- 
ness enhancements from ~500 to 800 au (Kalas et al. 2000; 
hereinafter KLSS), that have been attributed to a series of 
rings within the disc's midplane. Similar features in the SW 
extension were not detected and so eccentric rings were hy- 
pothesised. 

KLSS presented a dynamical model in which the stel- 
lar fly-by scenario (Kalas & Jewitt 1995) for generating the 
disc asymmetries was applied to the /3 Pic system. Their 
finding was that the asymmetry types present in the j3 Pic 
disc have analogues in the dynamical response of a particle 
disc to a perturbing stellar fly-by encounter. Furthermore, 
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they demonstrated that these features could occur simulta- 
neously with the formation of eccentric circumstellar rings. 
A simple analysis of the perturbed disc model was found to 
be in good quantitative agreement with the observations of 
both P Pic's disc asymmetries and the midplane features; for 
particular values of the viewing angle and epoch, and fly-by 
parameters. In this paper we present new observations and 
measurements of the length asymmetry of j3 Pic's scattered- 
light disc, and then go on to simulate the dynamics of a 
stellar flyby encounter with an initially symmetrical disc of 
particles. 

Similar numerical modelling to that presented here and 
in KLSS has been considered previously in the context of 
angular momentum transfer between the perturber's orbit 
and a protostellar disc (Clarke & Pringle 1993, Hall, Clarke 
& Pringle 1996), the dynamics of the Edgeworth-Kuiper belt 
(Ida, Larwood & Burkert 2000), and the tidal interaction 
of spiral galaxies (Toomre & Toomre 1972). Here we shall 
consider a wider range of encounter velocities (from bound 
orbits through to hyperbolic trajectories) and stellar masses 
(from M dwarfs to A stars, when scaled to the mass of f3 Pic) 
than has been studied before in any context. We also include 
a discussion of the orbital dynamics of captured particles 
with respect to the perturber. This paper can be considered 
as a companion to KLSS since here we also supplement their 
finding of circumstellar ring formation in similar simulations 
with an explanation of the dynamical origin of the features. 

As in other works (Whitmire, Matese & Tomley 1988, 
Mouillet et al. 1997, KLSS) we assume that the simula- 
tion particles' distribution represents the distribution of an 
underlying disc of planetesimals, which are parent bodies 
whose infrequent and therefore dynamically inconsequential 
collisions supply the dust scattering surface that is observed 
in the real system. A more complete model should include 
consideration of grain removal and generation processes in 
determining the appearence of the perturbed disc (reviewed 
by Backman & Paresce 1993). However, we defer treatment 
of those issues to future work, and focus here on the first- 
stage problem of the dynamics. This stage is important in 
deducing the mass and orbital parameters of the postulated 
stellar perturber, for which it is possible to perform star cat- 
alogue searches (e.g. Kalas, Deltorn & Larwood 2000), and 
we proceed to do this by examining the length asymmetry, 
and other measures of the disc response, as a function of stel- 
lar encounter parameters. However, we choose to examine 
the length asymmetry as the principal diagnostic measure 
of the disc response since it is the most distinct such feature 
both in the observations and the dynamical models. In our 
future work we shall also consider the subtler asymmetries 
in greater detail. 

In Section 2 we present deep optical coronagraphic im- 
ages of the disc. In Section 3 we describe the results of our 
numerical investigations of the dynamics of tidally disrupted 
particle discs, and provide general discussion of the induced 
asymmetries and particle stripping. In Section 4 we discuss 
our findings and summarise the conclusions. 



2 NEW OBSERVATIONAL RESULTS 

New optical coronagraphic observations of the environs of 
/3 Pic were obtained with the University of Hawaii 2.2 m 




Figure 1. PSF-subtracted, i?-band coronagraphic image cen- 
tred on f3 Pic. North is up and east is left. Bars in the north- 
east and southwest quadrants represent 10 arcsec and are placed 
103.6 arcsec (2000 au) radius from /3 Pic. The NE disc midplane is 
detected near 2000 au radius, whereas the SW midplane is traced 
to approximately 1000 au radius. Other large scale extended fea- 
tures in this image are due to instrumental scattered light (e.g. 
radial linear features, a narrow ring segment in the southeast 
quadrant, and amorphous structure near a bright star north of /3 
Pic). 

telescope on 2000 January 31. A key goal was to maximize 
sensitivity to the outer regions of the disc. The coronagraph 
used has 12 interchangeable, circular, focal-plane occulting 
masks, from which we selected a mask 12.5 arcsec in diam- 
eter. This particular mask size is approximately two times 
larger than the minimum size allowable by the seeing condi- 
tions, and it partially obscures the image out to 120 au pro- 
jected radius about /3 Pic. However, it also minimizes the 
scattered light in the field that dominates the background 
noise. Additionally, the occulting mask increases the effi- 
ciency of the observations by maximizing the integration 
time permitted before CCD saturation occurs just beyond 
its edge. 

We obtained 122 integrations of 30 seconds each 
through an i?-band filter (A ~0.67^tm). Re-imaging optics 
provided a pixel scale of 0.407 arcsec pixel -1 on a 2048 x 2048 
CCD. The field was masked by the coronagraphic focal plane 
assembly which has a circular clear aperture with diame- 
ter 336 arcsec. The airmass of observation (>3.0) combined 
with high winds, produced image quality varying from 1.0 
to 1.8 arcsec. For final data analysis, we selected 77 frames 
with the image quality characterized by stellar full-width at 
half-maximum (FWHM) < 1.2 arcsec. 

Data were bias-subtracted, flat-fielded, and sky- 
subtracted using the median sky value in the corners of the 
field for each frame. Comparison stars were also observed 
for the purpose of subtracting the /3 Pic PSF. However, this 
subtraction step adds noise and so we chose to remove the 
background using an artificial point-spread function (PSF). 
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Figure 2. Surface brightness contours for the NE midplane (top) 
and SW midplane (bottom). The original image (Fig. 1) has been 
rotated 60° anticlockwise such that the midplane lies along pixel 
rows. This image is separated along the vertical axis and the 
NE midplane then reflected such that the two sections of mid- 
plane point in the same direction. The outer contour represents 
25.0 mag arcsec -2 and has radial extent approximately 25 per 
cent greater for the NE extension than for the SW extension. 
The binning procedure for the data is described in the text. 



The artificial PSF is a figure of rotation for a seventh-order 
polynomial that is a least squares fit to j3 Pic's PSF sampled 
in radial directions perpendicular to the disc midplane. 

Fig. 1 shows the final image of /3 Pic with effective inte- 
gration time 2310 seconds. The 3er background noise corre- 
sponds to a sensitivity of R — +25.8 mag. The midplane of 
the SW extension is difficult to trace visually as a distinct 
linear feature beyond ~55 arcsec (1062 au), whereas the NE 
midplane is traceable at up to ~95 arcsec (1835 au) radius 
from 13 Pic. By these measurements the length asymmetry 
is approximately 70 per cent. 

A key problem in comparing the lengths of the two ex- 
tensions is that the morphology of the SW extension departs 
significantly from a linear feature, becoming vertically dis- 
tended and amorphous at large radii. To compensate for this 
difference in vertical structure we rotate the image such that 
the disk midplane lies along pixel rows and bin the data in 
the vertical direction. We do this by summing the light at 
each pixel location with the light from its two vertically- 
adjacent pixels. This vertical binning procedure smoothes 
the data and equalizes the NE and SW FWHM measure- 
ments when taking vertical cuts at given projected radii 
along the midplane extensions. Smoothed contours of the re- 
sulting image are shown in Fig. 2. Here, the 25.0 mag arcec -2 
contour traces the NE extension to 95 arcsec radius, and the 
SW disk extension is traced to ~75 arcsec (1450 au). The 
length asymmetry by this measurement is therefore ~25 per 
cent, which is fractionally larger than the length asymmetry 
measured by Kalas & Jewitt (1995) at about a half of the 
radial scale probed here. The asymmetric vertical morphol- 
ogy evident in the Fig. 2 contours also confirms the butter- 
fly asymmetry (see Section 3.3) reported by Kalas & Jewitt 
(1995). 



2.1 Summary 

The new data reported above show that the (3 Pic disk has 
a significantly greater radial extent (1835 au for the NE disk 
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Figure 3. The face-on and edge-on views of a sample of simula- 
tion particles in the computational frame centered on the primary 
mass (labelled Fig. 2 in KLSS). The primary is of unit mass and 
G = 1, giving a time unit such that the circular orbit period at 
unit radius is 27r. The horizontal axes correspond to the line of 
nodes of the perturber's orbit. Parameters are: perturber mass 
of 0.3, parabolic pericentre distance of 2.6, pericentre inclination 
of 30°. The disc of ~10 6 particles is initially distributed between 
0.2 and 2 radius units. The initial vertical density is exponentially 
decreasing, with an outer height of 0.6, and the aspect ratio and 
the radial surface density follow radial power laws of index —1.5. 
A random sample of ~25 per cent of the particles outside the ini- 
tial radius are shown in the top frame, the bottom frame shows a 
thin slice of particles about the horizontal axis. The path of the 
perturber is indicated with a dashed line and it enters the frame 
from the right-hand side. The particles' positions and velocities 
are integrated with a second-order leapfrog scheme; the perturber 
orbit is prescribed as a function of time. The configuration shown 
is for a time of 160 units after pericentre passage. The perturber 
was introduced at a separation from the primary of 20 units. 



extension; 1450 au for the SW) than has been previously 
measured. Smith & Terrile (1987) reported the maximum 
extent of f3 Pic's disc at 1300 au (i.e 67 arcsec) for the NE 
midplane at optical wavelengths. However, this was known 
as a sensitivity-limited value rather than the true physical 
extent of the disc. Here too we are unable to define the outer 
disc extents above the level of the background noise. 

The existence of a disc miplane out to ~10 3 au radius, 
as well as the ~10 2 au vertical thickness of the SW extension 
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Figure 4. Particle positions at various times for parabolic fly-by encounter simulations that use different components of the tidal 
potential. From the top row of frames to the bottom the times given are —3, +3, +20, and +60. The perturber is plotted with a cross, 
the primary star (origin) is plotted with an asterisk. In these simulations we used e = 10~ 7 , R p = 0.1, and R„ = 0.05, and removed 
particles that went beyond radii greater than 10 from the primary. 



(see Fig. 4 in Kalas & Jewitt 1995), is unexpected given our 
current understanding of planetesimal formation as applied 
to the Solar System (reviewed by Weidenschilling & Cuzzi 
1993). Thus we must conclude that planetesimals were cre- 
ated at much smaller radii than they are presently inferred 
to exist in (5 Pic, and subsequently suffered redistribution. 
Moreover this redistribution must allow for the existence of 
a lopsided and vertically flared disc. Global redistribution 
of this kind is a general outcome of a close stellar fly-by 
encounter (KLSS). Below we shall present numerical inves- 
tigations of the dynamics involved. 



3 DYNAMICAL MODELS 

In the fly-by model presented in KLSS a transient system 
of circumstellar eccentric rings was found to form in the 
perturbed disc of particles, as shown in Fig. 3. Fly-by pa- 
rameters, and viewing angle and epoch, were chosen to give 
the most satisfactory fit possible to the positions of the mid- 
plane brightness peaks found in the observations, and to the 
large-scale disc asymmetries. The model presented in KLSS 
will serve here as a reference case with which we compare 
the model behaviour as we vary the parameters about its 
particular values. In that model the computational length 
unit scaled to ~270 au with the initial disc radius (2 units) 
and distance of pericentre (2.6 units) for the perturber being 
~540 au and ~700 au respectively. Since G and the primary 
mass were set to unity in the computations, the computa- 
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tional time unit scaled to ~530yr (i.e. the period at unit 
radius divided by 2tt). Their system of computational units 
is used here in similar numerical simulations, however we 
point out that the scalings mentioned above are only used 
for reference, in order to maintain generality in the presen- 
tation of our results. 

The new numerical calculations presented below em- 
ploy a fifth-order accurate Runge-Kutta-Fehlberg integrator 
with a variable time-step size. While not as accurate as the 
Bulirsch-Stoer method (e.g. as used by Hall et al. 1996), this 
method allows for the fast integration of much larger num- 
bers of particles (~10 4-5 ) with moderate accuracy. Time- 
step size control is implemented as described in Press et al. 
(1992), by specifying a maximum relative error e in the ve- 
locity magnitude. The value of e was taken to be as small as 
practicable. Unless stated otherwise we take e -10" 9 . Thus, 
in a typical run of ~10 5 time-steps, the net relative error 
in the energy and angular momentum of any given parti- 
cle is ~1CP 4 . In practice the net relative error for most of 
the particles will be substantially smaller than this since the 
time-step size is controlled by a minority of particles that 
fall close to either of the stellar masses. To prevent too large 
a reduction in the time-step size owing to this effect we re- 
move particles from the calculation if they approach either 
star to within a radius R p of the primary or a radius i? s of 
the secondary. Unless stated otherwise we take: R p = 0.03 
and Rs — 0.01. The former value was chosen as the mini- 
mum required to ensure that the time-step size would not 
be controlled by close approaches of particles to the pri- 
mary rather than by close approaches to the secondary, for 
all model parameters considered here. The secondary was in- 
tegrated as if it were a test particle in motion about a fixed 
primary of mass equal to the total binary mass. The gravi- 
tational forces on particles are then prescribed as a function 
of position at each time-level with the addition of a fictitious 
force to compensate for the non-inertial origin (see below). 
As mentioned above, the test particles are collisionless and 
their distribution is taken as a tracer of the observed dust 
presumed to be produced by infrequent collisions in the real 
system. 

The origin of time is taken as the instant of pericentre 
passage, hence computational time-levels may be negative 
or positive. Throughout we shall use 10 4 test particles ini- 
tialised in Keplerian circular coplanar (zero inclination) or- 
bits, placed at random in the radial interval 0.5-2 such that 
the surface density profile decreases according to a power 
law of radius to the three-halves. The perturber has inital 
position and velocity components appropriate to the chosen 
orbital inclination, eccentricity, and mass. The perturber is 
initialised at a distance of ten times the initial outer disc ra- 
dius (i.e. at 20 units). We shall consider a single pericentre 
distance of q — 2.6 units, thus the initial disc size is ~0.8g. 
The length and mass units, and therefore the time unit, used 
in our calculations can be scaled so as to apply to any sys- 
tem within the range of dimensionless parameters that we 
cover. Thus, given our choice of relative disc size, fixing q 
does not affect the generality of our models. Also, since the 
initial disc is axisymmetric, the longitude of periastron for 
the perturber's orbit is arbitrary. In all our particle position 
plots the perturber orbit is orientated as in Fig. 3, except 
for differences in orbital inclination. 
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Figure 5. Data for the model using only the indirect component 
of the tidal potential (see Fig. 4). The top frame is for data taken 
at a time of +60; the bottom frame corresponds to a time of +180. 
Plotted are acos(/) against asin(/) for each particle, where a is 
the semimajor axis and / is the true anomaly. 
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Figure 6. Schematic representation of the initial disc (cross- 
hatched circle), and the orbital paths of a perturber for 
a parabolic encounter (solid line) and hyperbolic encounters 
(dashed lines, from lowermost to uppermost, for e = 2, 5, and 
10 respectively). The dotted line corresponds to a bound orbit 
with e = 0.75 (see Section 3.5). As for the simulations the origin 
is assumed to be fixed to the primary. 



3.1 Circumstellar ring formation 

In order to identify the dynamical origin of the ring features 
simulated by KLSS we consider the system response to dif- 
ferent components of the tidal potential (an approach also 
used recently by Pfalzner & Kley 2000, for gas discs). The 
total potential is given by ^ = + ^s, where the point 
mass potential at position r due to a primary of mass M p , 
to which we have fixed the origin, is: 
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(i) 



The tidal potential, due to the secondary with mass M s and 
instantaneous position r s , is given by: 



* = — 



GM S 



GM s r ■ r a 



(2) 



The first term on the right-hand side of (|2|) is the direct 
potential due to the secondary star; the second term is the 
indirect potential due to the acceleration of the origin. In 
a stellar fly-by this acceleration corresponds to the reflex 
motion of the primary star due to a gravitational tug from 
the secondary star. 

In Fig. 4 we give particle position plots for fly-by sim- 
ulations in which: we use only the direct term in the tidal 
potential, use only the indirect term in the tidal potential, 
use the full tidal potential. In the case that employs only the 
direct tidal term the disc is pulled towards the perturber as 
it approaches pericentre. In the case using only the indi- 
rect tidal term the disc is repelled by comparison as the 
perturber approaches pericentre (in the inertial frame this 
effect corresponds to acceleration of the primary towards the 
secondary) . In these first two cases eccentric nested rings are 
observed to form, although they appear as stronger features 
in the case where only the indirect tidal term is used. Since 
the effectiveness of the direct tidal interaction depends on 
the disc being able to closely approach the perturber, the 
role of the indirect tidal term when the full tidal potential is 
used is towards weakening the ring features that result from 
the direct tidal interaction, and to cause approximate peri- 
centric ant i- alignment of the two sets of rings. This scheme in 
the alignment of the particle orbits gives a length asymme- 
try between opposing sides of the disc when viewed edge-on, 
with the strongest ring features along the longest extension. 

In the above the features that appear to be ring-like are 
actually formed from distended tightly-wound one-armed 
kinematic spiral patterns. Different levels of advancement 
in the winding process can be seen at the third and fourth 
time-levels depicted in Fig. 4. This is more clearly illustrated 
in Fig. 5 where we plot particle positions according to the 
true anomaly (i.e. the angle between the instantaneous po- 
sition vector and the line of pericentre) for the last frame 
given in Fig. 4 and for a time much later after pericentre 
passage, for the model using only the indirect term in the 
tidal potential. [Although we have now identified the origin 
of the ring features, since over much of the azimuth the ap- 
pearence of the features resembles nested eccentric rings, we 
shall continue with our original usage of that description.] 

Winding of the kinematic spiral occurs as particles move 
around their orbits with local orbital frequency fl. Thus the 
time-scale for the generation of a ring after pericentre pas- 
sage is The lack of significant width to the spiral in 
Fig. 5 as compared to Fig. 4, indicates that the width of 
the rings is mainly due to the spread in the particle eccen- 
tricities. The rings will appear to propagate outwards and 
dissolve due to the growth of incoherence in the phase re- 
lationships of particles. Consider two neighbouring particles 
at radii or semimajor axes a a and ai, where a > ai. The 
decoherence time-scale id for the particle pair is given by: 



n(oi) 

2tt 



1 



(do) 



3/2 



(3) 



Table 1. Data for some fly-by encounter mod- 
els. The mass of the secondary relative to the 
primary is fi, the orbital eccentricity of the sec- 
ondary is e, and the orbital inclination of the 
secondary is i. The length ratio of extension A 
to extension B is L. We also give the percentage 
of the initial particle number that is ultimately: 
bound to the secondary star (n]-,), deleted by en- 
tering the secondary star's accretion radius (n s ), 
deleted by entering the primary star's accretion 
radius (n p ), or unbound from the system (n u ). 
The model listed in the third row was deemed to 
be most appropriate for modelling the asymme- 
tries in the Pic disc (KLSS, also Fig. 3). 
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In (5 Pic, and the simulations, the width of the midplane 
brightness features is ~5-10 per cent of their radii (KLSS). 
For dissolution of the ring features we require the tight wind- 
ing of a spiral down to the scale of a ring width, thus for a 
ring at radius r: 



td ~ 10 5 ( 

V500au/ 



3/2 



yr. 



(4) 



This agrees with the simulations of KLSS in which it is found 
that the ring features dissolve in about ten local orbital pe- 
riods. KLSS also found that viewing the perturbed model 
disc at approximately 10 5 yr after pericentre passage gave 
the most satisfactory fit to the observed disc, which is con- 
sistent with the innermost ring feature being detected at 
-500 au. 



3.2 Length asymmetry in coplanar models 

As discussed above, a length asymmetry in the disc is read- 
ily obtained from a close fly-by encounter. The length asym- 
metry in j3 Pic's disc is also the most easily measured of its 
asymmetric features. Below we shall continue to numerically 
investigate the form and size of the length asymmetry as a 
function of fly-by parameters. 



3.2.1 Parabolic versus hyperbolic trajectories 

A secondary star of orbital eccentricity e with respect to the 
primary star has relative velocity magnitude V q at pericentre 
distance q given by: 



l + e 



G(M P + M s 



(5) 



Hence for the scaling used by KLSS, eccentricities of 1, 2, 5, 
and 10 correspond to pericentre velocities of ~2.5, 3.1, 4.3, 
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Figure 7. The final distance of apocentre against the longitude 
of pericentre for fly-by simulations with perturber orbits of ec- 
centricity e = 1, 2, 5, and 10. Here and in similar plots we select 
particles with e < 1 and with separations from the secondary of 
greater than 10, at the final time-level. In all plots involving the 
longitude of pericentre, zero longitude is defined by the positive 
horizontal axis in face-on particle position plots. 



and 5.9 kms -1 respectively. In this Section we shall consider 
coplanar encounters of uniform pericentre distance 2.6, a 
perturber mass of 0.3, and these orbital eccentricities. The 
corresponding trajectories and the initial disc size used in 
these simulations are shown in Fig. 6. 

A useful diagnostic measure of the length asymmetry 
induced in the disc is given by plotting the final apocentre 
distance against the longitude of pericentre for each per- 
turbed particle orbit. By final we shall mean that the sec- 
ondary mass has passed through pericentre and has returned 
to a sufficiently large radius that its affect on the positions 
and velocities of the particles bound to the primary is negli- 
gible. In our simulations we integrated the system for a to- 
tal time of 120 units, which means that the final perturber 
distance was always much greater than the initial distance, 



-2 



longitude of pericentre 



Figure 8. The final distance of apocentre against the longitude 
of pericentre for a fly-by simulation with a parabolic perturber 
orbit in which only the direct tidal potential has been used. Data 
corresponds to those given in Fig. 4. 



corresponding to a time of at least 80 units after pericentre 
passage. This total run-time was also generally required in 
order for the numbers of particles either bound or unbound 
to each of the two stellar masses to stabilise (in practice this 
meant to approach a constant value with less than about 10 
per cent fluctuations; cf. Hall et al. 1996). This information 
is summarised in Table 1, for all the e > 1 models presented 
in this paper. 

In Fig. 7 we compare the disc response for encounters 
of different eccentricity. In all cases there are two oppositely 
directed extensions. The extension near to a longitude of 
—2 contains the strong ring features, and the other exten- 
sion (near longitude +1) contains the weak ring features, 
as discussed above. Respectively we shall refer to these as 
extension A and extension B. Thus extension A forms as 
a response to the indirect component of the tidal potential, 
and extension B forms as a response to the direct component 
of the tidal potential. 

As the encounter velocity (equivalently, the orbital ec- 
centricity) increases the peak of extension A diminishes in 
amplitude, while the peak of extension B increases its am- 
plitude. Further examination of the extension B results illu- 
minates a key difference between the parabolic (e = 1) and 
hyperbolic (e > 1) orbit cases. In the former, extension B 
is compact, but in the latter cases all show a much longer 
extension B. In Fig. 8 we show the final apocentre distance 
against the longitude of pericentre for the parabolic case in 
which only the direct part of the tidal potential is used. As 
we might have expected extension A is not in evidence (being 
due to the indirect part of the tidal potential), and exten- 
sion B is similar in nature to that found in the hyperbolic 
orbit cases. Thus we conclude that in the full tidal potential 
case the indirect component acts on the disc prior to peri- 
centre passage such that the disc response due to the direct 
component near pericentre is modified as compared with a 
calculation that uses only the direct term (also discussed in 
Section 3.1). Hence the more impulsive (i.e. hyperbolic) en- 
counters yield a response that is similar to that found for the 
parabolic case that uses the direct term only, being mani- 
fest in the form of extension B. It follows that higher velocity 
(i.e. more eccentric) perturbers will have extension B forms 
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Figure 9. The final eccentricity against semimajor axis distri- 
bution for fly-by simulations with perturber orbits of eccentricity 
e = 1, 2, 5, and 10. 
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Figure 10. The final apocentre distance against longitude of peri- 
centre distribution for parabolic fly-by simulations with perturber 
masses of 0.1, 0.3, 0.5, and 1.0. 



more like that found in Fig. 8, and smaller amplitudes for 
extension A. This can be seen in Fig. 7. 

In Fig. 9 we give the final eccentricities of the disc par- 
ticles for different eccentricities. In the parabolic orbit case, 
extension A particles correspond to the wing-shaped distri- 
bution having eccentricities approximately in the range 0.2- 
0.6. Superposed with this are the extension B particles that 
have generally very high eccentricites. The change in the 
relative lengths of extensions A and B with increasing per- 
turber eccentricity is also apparant. Generally, the excited 
eccentricities of the disc particles are lower for larger per- 
turber eccentricities, which can be explained from the sim- 
plistic impulse model for tidal interaction (Lin & Papaloizou 
1979): higher relative velocity results in weaker coupling be- 
tween the particles and the perturber. 



3.2.2 Variation of the perturber mass 

In Fig. 10 we show the final distance of apocentre against the 
final longitude of pericentre for parabolic encounters with 
different masses for the perturber. We consider the mass ra- 
tio fi = M s /Mp. Increasing the perturber mass from fj, = 0.1 
to (j, = 1.0, we find that extension A increases in ampli- 
tude and extension B diminishes in amplitude. The former 
effect is due to the larger perturber mass inducing a larger 
amplitude reflex motion of the primary star; the latter effect 
evidently results from the more massive perturber being able 
to strip more particles from the disc in the direct interaction 
that sets up extension B (see Section 3.4). 

The lowest mass model with fx — 0.1 has an extension 
B that has a similar extended form to the high eccentricity 
cases discussed above. This is due to the low mass of the 
perturber giving only a small reflex motion of the primary, 
and so a relatively weak interaction prior to pericentre pas- 
sage, that necessarily results in a more impulsive response. 
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Figure 11. The final distance of apocentre against the longitude 
of pericentre for parabolic fly-by simulations with orbital inclina- 
tion i = 30° and 60° . 

Only the /z = 0.3 model gives a compact extension B and 
extension A. This case has in terms of a model for the length 
asymmetry in the disc of /3 Pic the most appropriate am- 
plitude ratio for the extensions A and B of those coplanar 
models tried here. The ratio of lengths of apocentre for this 
model is approximately 9/6 = 1.5, the length ratios of the 
other models clearly being very different from unity. The 
length ratio between extensions A and B is given in Table 
1, for all our models with e > 1. Where the ratio is clearly 
much different from unity we simply state either > 1 or 
<C 1, since the interpretation of the length of an extension 
as defined by a small number of particles at large distances 
is not clear. 

We note that the size of extension A can be reduced by 
increasing the perturber's eccentricity, so it may be possible 
that a more massive (/x > 0.3) and eccentric orbit (e J; 
1) perturber could fit the observational results. However, 
since eccentricities larger than approximately unity do not 
display a compact extension B, we favour the perturber mass 
originally considered by KLSS (i.e. fi ~ 0.3, for e ~ 1). 

3.3 Inclined orbits 

In general, if the perturber's orbit and the initial disc plane 
are not coplanar, then the tidal interaction occurring be- 
tween the disc particles and the perturber can excite verti- 
cal motions in the particle orbits. Here we shall consider this 
effect in terms of the inclination i of the perturber's orbit to 
the initial disc plane. This angle is generated by rotating the 
perturber's orbit about the latus rectum (i.e. the line per- 
pendicular to the line of pericentre, that also passes through 
the origin), which then becomes the line of nodes (e.g. see 
Fig. 3). 

KLSS found that (for i — 30°) extension B became ver- 
tically flared, with much weaker inclinations in extension A, 
as can be seen in Fig. 3. From the preceeding discussion 
we deduce that the flaring results from the stronger direct 



tidal interaction, and the smaller inclinations generated in 
extension A result from the weaker indirect tidal interaction. 
The resulting morphology resembles the butterfly asymme- 
try (Kalas & Jewitt 1995, KLSS) in which the vertical height 
of the imaged disc is greatest north of the midplane in the 
southwest extension (see Figs. 1 and 2), and south of the 
midplane in the northeast extension, with the overall great- 
est vertical height in the former. 

3.3.1 Length asymmetry 

In Fig. 11 we give the final distance of apocentre against 
the longitude of pericentre for parabolic encounters with 
H = 0.3, and i — 30° and i = 60°. In the lower inclination 
case the disc response is very similar to the coplanar en- 
counter presented above, except that extension A is shorter 
by about 20 per cent. In the higher inclination case the re- 
sponse resembles that of a coplanar hyperbolic encounter, 
having a very short extension A and a very long extension 
B. As before, this is consistent with the notion that the 
disc interaction prior to pericentre passage is a key differ- 
ence between hyperbolic and parabolic encounters. In the 
i — 60° case the high inclination takes the perturber or- 
bit further from the disc, making the pericentre interaction 
more impulsive. This shows that we require an orbit with not 
too high an inclination, if we wish to preserve the compact 
nature of extension B. Additionally, the other angle possible, 
which corresponds to rotation of the orbit about the line of 
pericentre, need not be considered since for low values of i 
that angle takes the orbit yet further from the disc. In other 
words, we only consider an argument of pericentre of 90°. 

In the higher inclination case the length ratio of exten- 
sion A to extension B is clearly much less than unity. In 
the lower inclination case extension B is compact, but the 
peak is less distinct than for the coplanar case. We measure 
its height by taking the intersection of the tangents to the 
sides of the envelope of the extension particles near its base. 
This yields a height of approximately 6, which is similar 
to the coplanar model. Thus we deduce a ratio 7/6 ~ 1.2. 
This is the best fitting value, in terms of the length asym- 
metry, of those models presented in this paper. Scaling ac- 
cording to the KLSS model the lengths of the extensions 
are 1890 au and 1620 au, which is in good agreement with 
the new determination of the size of the disc presented in 
Section 2. 

The exact numerical value of the length ratio deter- 
mined for the P Pic disc could be recovered in our mod- 
els by taking a slightly smaller orbital inclination than 30° . 
That would not be productive at this stage, however, since 
here we are not concerned with the details of translating 
our results into a scattered light disc model, which could in- 
troduce a small correction factor. However, for reference, a 
model with i — 15°was tried and was found to give a length 
ratio of 1.4. Thus we expect the appropriate inclination to 
lie between 15° and 30°, i.e. closer to the lower inclination 
model above than to coplanar. 

3.3.2 Flaring and warps 

In Fig. 12 we show the final eccentricity and inclination dis- 
tributions for the disc particles in the two cases mentioned 
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Figure 12. The final eccentricity and inclination against semimajor axis distributions for parabolic fly-by simulations with orbital 
inclination i = 30° and 60° . 




Figure 14. The final eccentricity and distance of pericentre against semimajor axis for a hyperbolic fly-by simulation with orbital 
eccentricity e = 2 and orbital inclination i = 30°. 



above. The amplitude of the response is lower for higher 
values of i, for the reasons indicated above. The eccentric- 
ity response resembles that of a coplanar hyperbolic orbit, 
as expected from the preceeding discussion. The inclinations 
associated with extension A are much smaller than for those 
associated with extension B. However, the extension A in- 
clinations are higher in the i — 60° model, owing to input 
from the direct interaction because of increased impulsive- 
ness of the interaction occuring at the ascending node in 
that case. Thus, by symmetry of the orbit with respect to 
the initial disc, we expect the inclinations in each extension 
to become comparable when orthogonality is approached by 
increasing i. This provides a further need for a low inclina- 
tion encounter as the disc of /3 Pic is strongly lopsided with 
respect to vertical flaring (Kalas & Jewitt 1995). 

For these two cases we also give the final pericentre dis- 
tance distributions in Fig. 13. In the low inclination case we 
see that particles from extension B can plunge to pericentre 
distances of 0.2. KLSS proposed that these particles could 



explain the existence of an inner 'warp' in the scattered light 
dust disc of P Pic at ~50 au (Burrows et al. 1995, Mouillet 
et al. 1997). We can see from Fig. 13, that for the higher in- 
clination encounter, extension B particles have significantly 
higher perturbed pericentre distances (down to about 0.7; or 
~200 au). This is a direct consequence of the lower eccentric- 
ities produced at comparable semimajor axes in extension B 
for the two cases. Thus the minimum pericentre distance is 
minimised in coplanar cases. Correspondingly, larger mini- 
mum pericentre distances are also charateristic of hyperbolic 
encounters. Fig. 14 gives the final eccentricity and pericen- 
tre distributions for a hyperbolic encounter with e = 2 and 
i = 30°. If the explanation of the inner warp suggested by 
KLSS is confirmed then this also supports the need for a low 
inclination and relative velocity encounter. 
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Figure 13. The final distance of pericentre against semimajor 
axis for parabolic fly-by simulations with orbital inclination i = 
30° and 60°. 



3.4 Planetesimal capture 

As can be seen in Fig. 4 the close fly-by of a disc of particles 
by a perturbing star can result in the very close approach of 
some of the disc particles to the perturber. Since the unper- 
turbed disc does not extend all the way out to pericentre, 
this occurs by particles streaming through the vicinity of 
the instantaneous Lagrange Point between the two stars, 
about the time of closest approach. In the case that parti- 
cles become captured by the perturber the capture efficiency 
is generally low (see also Hall et al. 1996). Table 1 gives the 
percentage of initial disc particles that are captured by the 
perturber for the models discussed above, and from it we can 
see that the lower eccentricity and inclination encounters are 
required, in order to effect capture at all. 

For the model parameters used by KLSS, the percentage 
of captured particles is ~10 per cent of the initial particle 
number, but note that if the initial disc were allowed to ex- 
tend to pericentre this figure could be fractionally larger 
(Hall et al. 1996). In the case of /? Pic, considering the 
maximal parent body disc of Backman & Paresce (1993), 
the captured material could therefore represent up to about 
10 3 earth masses of planetesimals which would subsequently 
orbit the perturber and produce a scattered light signa- 
ture similar to Vega-like sources. This would not however 
uniquely identify a perturber candidate since approximately 
15 per cent of main sequence stars of all types are known to 
have such discs (Backman & Gillett 1987, Aumann 1988), 
but given the encounter parameters we find here, it seems 
likely that a perturber candidate should have this signature. 

None the less we may find crude approximations for the 
orbital properties of the captured planetesimals. We shall 
make an analogy with the theory of mass transfer in semi- 
detached binaries, in which the captured matter orbits the 
secondary with the specific angular momentum of relative 
motion at the Lagrange Point (see Frank, King & Raine 
1992). Here we deduce that if a particle in the outermost 



parts of the disc is captured close to the perturber's peri- 
centre, then it will orbit the perturber with specific angular 
momentum: 



G(l + e)(M p + M s ) 



GM P 
Ro 



(6) 



where R a is the initial outer radius of the disc and b is the 
separation of the secondary from the Lagrange Point at the 
instant of pericentre passage, approximated by: 



- = i +0.227 log( M ). 
q A 



(7) 



In the above we have taken the circular velocity of a particle 
lying on the binary axis relative to the perturber at pericen- 
tre for the input particle velocity. Although the instanta- 
neous Lagrange Point is located well within the envelope of 
the initial disc, the rapid eccentricity changes in the disc 
prior to particle transfer cause compression of the disc such 
that the outer disc particles pass through the vicinity of the 
Lagrange Point with approximately their initial velocities 
(e.g. see Fig. 4). 

With respect to the perturber, the specific angular mo- 
mentum of the captured particle is: yj GM a a'(l — e' 2 ), where 
the characteristic semimajor axis and eccentricity with re- 
spect to the secondary are respectively a' and e'. Hence we 
have: 



a(l-e )=g(l + e) 1 



1 - 



q/Ro 



(l + e)(l + /i) 



(8) 



Thus for a model with e = 1, fj, = 0.3, and q = 2.6, we find 
a' (I — e' 2 ) ~ 0.3. This equation does not apply when the 
relative velocity of the perturber and disc particles is suffi- 
ciently high that particle transfer is energetically disallowed. 
As shown in Table 1 we find that particles are captured when 
e < 2 but not when e > 5. 

The above is in fair agreement with the simulations. In 
Fig. 15 we give the positions, eccentricities, and semimajor 
axes of captured particles relative to the perturber for the 
above model. We find the bound particles have a median 
value a' ~ 0.70, from which we calculate e' ~ 0.76. The 
actual median value is e' ~ 0.85 for the particles, which is 
~ 10 per cent larger than predicted. This is likely due to the 
continued interaction of the primary and the captured par- 
ticles as the perturber receeds from pericentre. The particle 
transfer does not occur in a well defined stream which re- 
sults in a large spread in semimajor axes and eccentricities. 
However, the captured particles appear to manifest as an 
asymmetrical disc about the perturber (upper panel of Fig. 
15), and show clustering in e' and a' (lower panel of Fig. 15), 
indicating some level of coherence in the transfer (seen in 
Fig. 4). Scaling according to KLSS, we find that 90% of the 
disc particles are located within 500 au of the perturber. For 
the KLSS model (i.e. the previous model augmented with 
an orbital inclination of 30° for the perturber), the median 
eccentricity and semimajor axis are each found to be ~ 0.7, 
consistent with the higher relative disc velocity expected for 
the inclined encounter. 

The luminosity of an MOV dwarf star (of mass ~ 



© 0000 RAS, MNRAS 000, 



12 J.D. Larwood and P.G. Kalas 




Figure 15. The final distance of apocentre against the longitude of pericentre for a bound orbit fly-by simulation with e = 0.8. The 
companion starts at apocentre and was set up and integrated much as for the hyperbolic orbit cases. Data is given for the first and 
second complete orbits. In this simulation we used e = 10~ 7 . 



0.5 Mq, corresponding to fi = 0.3) is ~10 2 times smaller 
than for f3 Pic. If the amount of scattering surface in the 
disc is proportional to the square of the number of parent 
bodies (e.g. Backman & Paresce 1993), then as an upper 
limit the brightness of the perturber's disc is expected to be 
~10 4 times fainter than that around fJ Pic, and ~10 2 times 
fainter relative to the central star. Thus it is plausible that 
the perturber may demonstrate a measurable scattered light 
disc as a result of collisional erosion amongst the captured 
planetesimals at the distance to j3 Pic. 




3.5 Close encounters with a bound perturber 

Since the innermost midplane brightness feature identified 
by KLSS is measured to be at ~500au the age of the rings 
in the fi Pic disc should be ~10 5 yr, as calculated above. 
Although this figure could be revised up or down according 
to future observations, we shall take this to mean that the 
epoch of the fly-by encounter was about 10 5 yr in the past. 
Use of the Hipparcos survey data reveals that this is also the 
minimum time-scale we must look back from the present in 
order to find a perturber fly-by at ^ 10 au (Kalas et al. 
2000). The eighteen most closely approaching stars identi- 
fied by Kalas et al. (2000) have a mean encounter velocity 
~40kms _1 , which is considerably higher than that deduced 
from our dynamical calculations. Thus we are prompted to 
consider companions that are or were physically associated 
with fj Pic. 

Repeated close encounters of the bound star with the 
disc may destroy the asymmetrical disc structure and ring 
features set up by the first encounter. Since for e < 1 the 
semimajor axis of the perturber can be written q/(l — e), by 
demanding the orbital period to be greater than some value 
fmin, we deduce: 



e > 1 



4tT 



G(M P + M S )P2 



1/3 



(9) 



Thus if we take P m i n = 10 yr we infer that the perturbed 
binary eccentricity needs to be greater than 0.75, in order 
that we may view the ring features at their inferred age. 



Figure 16. Particle positions and eccentricity versus semimajor 
axis for captured particles relative to the perturber. 



In Fig. 6 we show a trajectory with e = 0.75. The simi- 
larity of such high eccentricty bound orbits to the parabolic 
trajectory suggests that the disc response should be similar. 
In Fig. 16 we give the apocentre distance versus longitude 
of pericentre plots for a coplanar encounter simulation with 
e = 0.8, fi — 0.3, and q = 2.6. The perturber starts at apoc- 
entre and we give results for the first and second completed 
orbital periods. Clearly the length asymmetry is completely 
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disrupted after the second encounter. In Fig. 17 we give the 
corresponding particle positions, which show that the ring 
structure is also destroyed after the second encounter. Thus 
as suspected only a single fly-by can be allowed. 

Possible scenarios in which a single pericentre passage 
could have occurred in the recent past are: 1) the perturber 
was originally a bound companion to /3 Pic with a sufficiently 
large pericentre distance that no significant interaction with 
the disc occurred over its lifetime, and sufficently large semi- 
major axis and eccentricity that it could be perturbed into 
the required trajectory by the close fly-by of a massive field 
star ~10 5 yr in the past; 2) (3 Pic was a hierarchical multiple 
system that was dynamically unstable on a time-scale com- 
parable to its present age, resulting in the ejection of one or 
more stars and the development of a disc-intercepting orbit 
for another, or an ejected star passed by /3 Pic on its way 
out of the system. 

In the first scenario the three-body interaction should 
have occurred at least 10 5 yr ago. The further in the past we 
put this event, the more eccentric the resultant binary orbit 
must be in order to delay closest approach until 10 5 yr ago (if 
the perturbed star has an initially out-going trajectory). Al- 
though perhaps more attractive than the lower probability 
close binary encounter model used by KLSS, the problem of 
low relative velocity is shifted from the star-disc interaction 
to the interaction between the second and third bodies. In 
respect of the second scenario we point out that Weinberger 
et al. (2000) report an age of 5 Myr for the triple system 
HD 141569. This consists of a Vega-like BVe primary with 
a 500 au disc and two M stars each at ~1000au projected 
separation. Weinberger et al. (2000) deduce that the rel- 
ative separations involved may make the system unstable, 
giving the possibility that HD 141569 represents an earlier 
stage in the evolution of the (3 Pic system (with an age 
of 20 Myr). Studies of the dynamics involved suggest that 
young stars with multiple companions may not lose some 
of them via dynamical instability until main sequence ages 
(Eggleton & Kiseleva 1995), which is supported by obser- 
vations of high binary frequency amongst young stars com- 
pared to field stars (Mathieu 1994). 



4 DISCUSSION 

Here we have considered the stellar fly-by hypothesis that 
may account for the asymmetrical structure in the scattered 
light dust disc of (3 Pic. We investigated the dynamics of 
the perturbed planetesimal disc in an attempt to learn more 
about the possible parameters of a fly-by encounter. Addi- 
tionally, we have described the origin of transient circum- 
stellar eccentric ring structures as being a general outcome 
of an encounter, due to the reflex motion of the primary 
star as the perturber passes through closest approach. The 
apparent ring system is actually an eccentric tightly-wound 
one-armed spiral pattern that gradually disappears owing 
to phase-mixing in the particles' orbital motion. In inclined- 
orbit simulations we find that the vertical disc response also 
has analogues in the scattered light disc of (3 Pic. We sum- 
marise the coplanar results in Fig. 18, and the inclined-orbit 
results in Fig. 19. 

For the perturber we conclude that the mass should 
not be very different from 0.5 Mq (corresponding to spec- 



first return to epastron secsnd return to apestron 




Figure 17. Particle positions corresponding to data given in 
Fig. 16. The primary's position is plotted with an asterisk. 



tral type MOV). In this case we favour a low inclination pro- 
grade encounter with relative velocity of a few kms -1 , which 
makes a bound perturber seem the most natural choice. One 
problem with this is that asymmetries would be destroyed 
in the second pericentre passage. Thus, given the deduced 
youth of the asymmetries, the perturber in this scenario is 
required to make a single close approach in the recent past 
(i.e. ~10 yr ago). This can be accounted for by making the 
perturber-/? Pic system an initial wide binary that is per- 
turbed by a massive passing field star ~10 5 yr in the past, 
such that the perturber is sent into a close approach orbit 
about (3 Pic. A more probable solution is to arrange for f3 Pic 
to have originally been a hierachical multiple system (cf. HD 
141569; Weinberger et al. 2000) that became dynamically 
unstable on the age of the star, and with the result being 
that one component was sent into a disc-intercepting orbit. 
In any case the perturber should have captured planetesi- 
mals from the primary's disc and may exhibit a detectable 
scattered light disc as a result. A survey of a more extended 
field around [3 Pic than has been considered previously is in 
progress in order to test our hypotheses. 

The dynamical response of the models suggests that 
the (3 Pic disc might be separable into different components 
corresponding to groupings of perturbed particle orbits. For 
example, we find that outwards of ^250 au the SW exten- 
sion may not have a distinct midplane owing to pumping of 
particle inclinations. The positions close to pericentre of the 
moderately eccentric NE midplane particles could account 
for the measurement of a SW midplane out to ~650 au, in 
superposition with the more diffuse vertically and radially 
extended component. When the optical image is vertically 
compressed (Fig. 2) we are able to measure the SW extension 
to much larger radii (1450 au). In our model the correspond- 
ing (very eccentric and inclined) particle orbits also intersect 
the disc near to their pericentres (down to ~50 au) possibly 
accounting for the inferred inner warp of the disc, which is 
aligned with the outer flared envelope in the SW. In sum- 
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Figure 19. Summary of inclined-orbit models in face-on and 
edge-on particle plots. The models were run for 120 time units, 
with the perturber initialised a radial distance of 20 from the pri- 
mary. The left-hand panel shows results for the model with e = 1, 
fi = 0.3, and i = 30°. The right-hand panel shows results for the 
corresponding case with i = 60°. Uppermost frames are face-on 
views, middle frames are edge-on projections along the horizontal 
axis (line of nodes), the bottom frames are edge-on projections 
along the axis defining i = 0. Thus the projections are on to the 
three mutually orthogonal Cartesian planes of the computational 
coordinate system. The position of the primary is plotted with an 
asterisk. 



Figure 18. Summary of coplanar model results in face-on particle 
plots. The models were run for 120 time units, with the perturber 
initialised a radial distance of 20 units from the primary. After 
much longer times the phase-mixing in the particles' orbits results 
in dissolution of the ring features, maintaining the distended en- 
velopes of the discs seen above. The position of the primary is 
plotted with an asterisk. The top left-hand side panel shows re- 
sults for models with /i = 0.3 and from top to bottom e = 2, 5, 
and 10 respectively. The top right-hand side panel shows results 
for models with e = 1 and from top to bottom /i = 0.1, 0.5, and 
1 respectively. The bottom panel shows results, on a larger scale, 
for the model with e = 1 and fi = 0.3. The perturber trajectory 
is shown with a solid line. We see: the perturbed disc with rings, 
the captured particles clustered about the instantaneous position 
of the perturber, and escaping particles (the large arcs on the 
left-hand side of the plot). 



mary, there are three recognizable particle groupings that 
can be related to the morphology of the (3 Pic disc. These 
are: the highly eccentric and inclined particles that reach 
apocentre in the SW (extension B particles) , the moderately 
eccentric and weakly inclined particles that reach apocentre 
in the NE (extension A particles), and the relatively unper- 
turbed particles inside ~200au radius (Figs. 12 and 13). 
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